Exotic Superfluid States of Lattice Fermions in Elongated Traps 
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We present real-space dynamical mean-field theory calculations for attractively interacting 
fermions in three-dimensional lattices with elongated traps. The critical polarization is found to 
be 0.8, regardless of the trap elongation. Below the critical polarization, we find unconventional su- 
perfluid structures where the polarized superfluid and Fulde-Ferrell-Larkin-Ovchinnikov-type states 
emerge across the entire core region. 



The nature of pairing in spin-polarized fermion sys- 
tems is a fundamental problem in many areas of 
physics, including superconductors in a strong mag- 
netic field, neutron-proton pairing in nuclear matter, 
and color superconductivity in high density QCD [D- 
3]. Non-BCS pairing mechanisms have been proposed 
for spin-polarized fermions. The Fulde-Ferrell-Larkin- 
Ovchinnikov (FFLO) states exhibit finite- momentum 
pairing, causing spatially oscillating pair potentials [J, Q . 
In the polarized superfluid states, so-called Sarma or 
breached-pair (BP) states, zero-momentum pairing oc- 
curs, having excess unpaired particles In trapped 
ultracold Fermi gases [ll|, [l2[ , recent experiments with 
spin imbalance observed the superfluidity and the tran- 
sition to the normal state with increasing polarization 
I3I- IB]. These systems offer an unprecedented access 
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to the search for exotic superfluid states. In this Letter, 
we present the first large-scale calculations beyond mean- 
field theory, explicitly considering the trap geometry. We 
characterize the critical polarization of the transition and 
show evidence for exotic superfluidity. 

The FFLO and polarized superfluid states have been 
originally proposed for translationally invariant systems. 
The FFLO state, in particular, has been suggested to 
be stabilized by reduced dimensionality. Therefore, it 
is important to incorporate broken translational symme- 
try caused by the presence of a trap in ultracold gases. 
Highly elongated traps were used in the experiments with 
spin-imbalanced 6 Li gases Concerning previous 

mean-field theory calculations with explicit trap confine- 
ment [TH , the FFLO-type oscillations have been expected 
to occur only at the narrow edges of the superfluid core 
in a spherical trap [19| . Quantum fluctuations are largely 
neglected in these mean-field calculations, and a full-scale 
ah initio approach remains elusive so far in three dimen- 
sions. Here we address this issue via an extension of dy- 
namical mean-field theory considering both the full local 
quantum fluctuations and the trap effects. 

We study the problem of a trapped lattice fermion sys- 
tem at zero temperature (T — 0) in three dimensions 
(3D) by usin g a real-space dynamical mean-field theory 
(DMFT) [2fj| - [23l |. With various trap aspect ratios ex- 
amined for elongated traps, we find the critical polariza- 
tion near 0.8 insensitive to the trap aspect ratios. In a 
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FIG. 1: Cross-sectional visualization of the particle density 
and pair potential structure in the elongated trap with the 
aspect ratio a = 7.5 for selected polarizations P < 0.8. (a) 
The density difference Sn = n-f — n± is plotted in the y = 
plane. The lighter (darker) gray indicates the larger (smaller) 
value of Sn. Black defines 8n = 0. (b) The pair potential 
A(x, y = 0, z) shows the evolution of the superfluid core struc- 
ture as a function of P. The FFLO-like oscillations are ob- 
served at high P's but below P = 0.8 where A vanishes. 



wide range of polarizations, we observe that polarization 
and finite pair potential coexist at the core region, indi- 
cating the polarized superfluid core, surrounded by the 
normal unpaired particles. In highly elongated traps, 
as we approach the critical polarization, we find that 
the FFLO-type states emerge with spatially oscillating 
pair potential across the trap center. For comparison, 
we present also Bogoliubov-de Gennes (BdG) mean-field 
calculations. 

Our real-space DMFT solves the attractive Hubbard 
model with an anisotropic trapping potential at T = 0, 

H = -t'y] c\ a c 3<y -u^2 n n n u + - av)««7, 

(ij)a i i<? 

where cl a (ci a ) creates (annihilates) a fermion with spin 
a at site i, the density operator m a = c\ a Ci a , and fi a 
denotes the chemical potential. The hopping t between 
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neighboring sites (ij) is set to unity. The onsite inter- 
action U is chosen as the unitarity value for the cubic 
lattice, U ~ 7.915 0. The trapping potential is given 
as Vi — Vo[x' 2 + y' 2 + (z'/a) 2 ] where Vq and a are the 
trap strength and aspect ratio. In the 3D lattice with 
size L x x L y x L z , the coordinates of sites are assigned as 
£' = £ - |, where £ 6 {x,y,z}, and £ = — . .., ^ - 1. 

The real-space DMFT is constructed for inhomoge- 
neous systems by considering local, yet site dependent, 
self-energy terms, which explicitly include trap effects be- 
yond the local density approximation (LDA). This ap- 
proximation in the DMFT allows s-wave pairing [20|. 
The impurity problem now becomes site-dependent and 
must be solved for each site with an efficient solver. We 
employ the exact diagonalization (ED) method to solve 
a 7-orbital impurity Hamiltonian with a superconducting 
bath [2(| which provides a good convergence (cf. (25[) in 
our reliability tests. There are two advantages in choos- 
ing the ED method: it is easy to consider the general- 
ized Anderson model required to describe the superfluid 
phases; it efficiently finds the ground state. The possible 
source of error is the effective bath of a limited number 
of orbitals, causing spiky spectral functions. 

We consider harmonic trapping potentials with various 
trap aspect ratios a = 1, 2.5, 5, 7.5, 10. The chemical po- 
tentials are adjusted to fix the total particle number N = 
A/ t + iV; ~ 210 with the polarization P = ( N t - N±)/N. 
The particle density is below quarter filling at all lattice 
sites. The trap strength is given as Vq — 0.02a 2 / 3 to keep 
the volume constant. The lattice size varies accordingly 
with a: (L x ,L y ,L z ) = (24,24,24) for a = 1, (16,16,40) 
for a = 2.5, (14, 14, 70) for a = 5, (14, 14, 80) for a = 7.5, 
and (14, 14, 100) for a = 10. The lattice size becomes 
a bottleneck in the real-space DMFT. The computation 
took ~ 10 9 CPU seconds in supercomputers. 

Figure [TJ shows how the distribution of particles and 
their superfluid characteristics change with polarization 
P in an elongated trap with aspect ratio a = 7.5. The 
density difference Srii = — and the pair poten- 
tial Aj = (cLct.) characterize the superfluid structure: 
for small P's, the particles at the trap center are fully 
paired superfluid (dark gray, 5n = 0, finite A) and sur- 
rounded by unpaired particles (light gray, 5n ^ 0). The 
small bumps in A are found along the axial {z) direction, 
which can be interpreted as an analog of the proximity 
effect at a superconductor-ferromagnet interface [26| . As 
P increases further, we find that two features become no- 
ticeable in A. For P > 0.7, an oscillatory structure in 
A emerges across the trap center, indicating that FFLO- 
type states may exist in this elongated system. When P 
reaches 0.8, a transition to the normal phase occurs, as 
indicated by the complete suppression of A. 

We study the transition to the normal phase system- 
atically with the order parameter, -q = J^. |Ai| 2 /7V^, av- 
eraged over the whole system. In Fig. Ufa) , it turns out 



(a) 



B 
<E 
E 

CO 



(b) 



o 
to 

N 

is 

o 

CL 



CD 

o 



0.5 
0.4 
0.3 
0.2 
0.1 


0.4 
0.3 
0.2 
0.1 





i o B soon 



cc=10.0 ■ 
a=7.5 
a=5.0 * 
a=2.5 
a=1.0 
LDA o 



X\Af I Nj, 




0.2 0.4 0.6 0.8 1 



0.2 



0.4 0.6 
Polarization (P) 



0.8 



FIG. 2: Transition to the normal phase with increas- 
ing polarization P examined for various trap aspect ratios 
a. (a) Superfluid order parameter ]T\ |Aj| 2 /iVi. vanishes 
above P c — 0.8, regardless of a. (b) Central polarization 
(n-f — ni)/(n^ + n{) of the site at the cloud center becomes 
finite far below P c , implying the presence of a polarized super- 
fluid core with finite 5n and A as indicated in the inset. The 
conventional DMFT data with LDA are given for comparison. 



that the behavior of f](P) and the critical polarization 
P c ~ 0.8 is insensitive to the trap aspect ratios that we 
have examined. Moreover, we find that our P c agrees 
with P c ~ 0.78 found in the experiments [l5| and [l6j |. 
and the quantum Monte Carlo (QMC) results in contin- 



uum with LDA 27J, l28( . The order parameter is compa- 
rable to the condensate fraction measured in 15]. 

On the other hand, we find that pairing in the super- 
fluid core in our systems is unconventional. Figure [2Kb) 
indicates the coexistence of finite polarization and order 
parameter at the center. This characterizes the polar- 
ized superfluid core and is apparent at intermediate or 
even at low P's for large a. In contrast to the order pa- 
rameter, the behavior of the central polarization depends 
on a. The LDA calculations, based on the conventional 
DMFT for homogeneous lattices at T = 0, show a first- 
order transition from the fully paired superfluid phase to 
the partially paired normal phase without the intermedi- 
ate polarized superfluid phase. However, the LDA can- 
not consider the inhomogeneous exotic phase observed in 
Fig. [TJ and it ignores the interface effects that may con- 
tribute to the polarized superfluidity. While the potential 
role of finite-size effects should not be ruled out, we have 
tested also N ~ 320 and found the same features. 

Let us discuss in detail how the cloud structure and 
the pairing evolve as P increases. Figure [3] presents the 
axial profiles of densities n^ t i, their difference 8n, and 
the pair potential A along the z-axis. The data with 
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FIG. 3: Evolution of the cloud structure with increasing polarization P. Axial density profiles n and pair potential A are 
shown along the z-axis. The plots are selected for the trap aspect ratio a = 7.5. The coexistence of the finite dip in the density 
difference Sn = — n± and a finite A at the central region characterizes the polarized superfluid core. The central dip in Sn 
disappears at high P — 0.74, replaced by the FFLO-type oscillations in A(z) spreading across the trap center. 



a = 7.5 are shown as a representative example. The 
fully polarized edges are well separated from the super- 
fluid core with the partially polarized intermediate re- 
gions, while the central region shows nonuniform A that 
changes with P. Our LDA data (not shown) gives the 
three-shell structure: the fully paired superfluid core, the 
partially polarized normal state shell, and the fully po- 
larized edges. This perfectly agrees with those from the 
QMC with LDA ^2j$. While, at low P's, our real-space 
DMFT presents a similar shell structure, the density pro- 
files indicate the exotic superfluid phase at high P's that 
is not accessible by the LDA with conventional DMFT. 

At very low polarization P = 0.08, the core is fully 
paired, and small oscillations of A appear in the partially 
polarized intermediate shoulders of Sn. As P increases, 
polarization becomes finite at center, leading to the wide 
dip area with finite Sn, indicating the polarized superfluid 
core. The dip region gets narrower with increasing P 
while the shoulders with small A oscillations get wider. 

When P increases further, the dip in Sn finally disap- 
pears, and the FFLO-type oscillations in A spread across 
the trap center, as plotted for P = 0.74 in Fig. [3l As P 
approaches P c , the amplitude of A becomes smaller and 
finally vanishes at P c . Along with the oscillating A, we 
also find the corresponding oscillations in Sn with a half 
period of A. This is consistent with the prediction of the 
FFLO-type phase in one dimensional (ID) systems [29| . 
However, our systems still show 3D features: at low P's, 
density profiles are similar to the LDA data in 3D while 
fully paired edges together with the FFLO phase at the 
center are expected in the strictly ID system 17, 29| . 

We find that these FFLO-type oscillations emerging 
across the trap center, rather than residing just at the 
edges, are found only in highly elongated traps with the 
trap aspect ratios a > 5.0 among the examined values of 
a's. The dependence of the core phase on a is summa- 
rized in Fig. [4] For a — 1.0 and 2.5, we have found the 
oscillations of A reside only at the edges for all P < P c . 

The diagram in Fig. @] separates the phases of the cen- 
tral region of the traps into fully paired superfluid (SF), 



polarized superfluid (pSF), FFLO, and normal (N). The 
FFLO-type phase can be signaled by the disappearance 
of the central dip in Sn, as shown in Fig. [3] The bound- 
ary between SF and pSF is rather unclear because Sn be- 
comes finite gradually with increasing P. However, our 
calculations show a tendency for the central polarization 
to grow faster as the trap gets more elongated. 

The nature of the polarized superfluidity found in our 
real-space DMFT calculations at T = can be very sub- 
tle. It has not been clear yet how stable pSF is at T = 0, 
especially in a trap. The BP state is in general ener- 
getically unstable toward the phase separation at a weak 
coupling at T = 3(1 31 1 . For a strong coupling, an ear- 
lier DMFT calculation in half-filled homogeneous cubic 
lattices found that the pSF may be stable down to a very 
low temperature (32[, while other calculations suggested 
that the density difference may decay exponentially with 
decreasing temperature in quarter-filled infinite dimen- 
sional lattices [33(. In our LDA calculations, we have 
observed a sharp SF-N transition without an intermedi- 
ate pSF phase. The stability issue can be more complex 
in our trapped asymmetric systems. We have found non- 



J3 
o 
CO 

ir 

o 

CD 
Q_ 
GO 
< 



10.0 

7.5 

5.0 

2.5 
1.0 



SF 



pSF 



FFLO N 



0.2 0.4 0.6 0.8 1 

Polarization (P) 

FIG. 4: Trap aspect ratio a and polarization P dependence 
of the phase of the central region of the trap. Fully paired 
superfluid (SF) is found at low P's and continuously evolves 
into the polarized superfluid (pSF). The pSF is specified with 
the central polarization larger than 0.005. At higher P's, the 
FFLO-type phase is identified for a > 5.0. The dashed lines 
between the phases are drawn for guidance. 
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FIG. 5: BdG calculations at unitarity in the elongated trap 
with q = 10. (a) The condensate fraction is plotted as a 
function of polarization P. (b) The density difference 5nedG 
and (c) the pair potential Ab<ig (normalized by the Fermi 
energy Ef) in the y — plane are shown for selected P's. 
The majority population Nf is fixed at 200 for a given P. 



standard (pSF)-FFLO-N type structures in the elongated 
traps, which cannot be explained with LDA. 

Our real-space DMFT inspires the question: how dif- 
ferent are the results qualitatively from those given by 
the BdG static mean-field theory, incorporating a trap- 
ping potential 34]? The DMFT is based on a discrete 
lattice system, nevertheless we choose to compare it to a 
BdG approach without a lattice (the density of our lattice 
fcrmions is well below the quarter filling) . In Fig. O we 
observe very high P c (~ 0.95) and the strong deformation 
of the core shape at high P's. These features are in clear 
contrast with those found in our real-space DMFT. The 
DMFT accurately describes local quantum fluctuations 
which are missing in the BdG theory that does not prop- 
erly include interactions in the normal state. Such effects 
cause the large deviation in P c in the strongly interacting 
regime. 

The explicit test for the trap aspect ratio sheds light 
on the comparison between the experiments HI, ll| and 
14j . Despite the difference between dilute gases and low- 
filling lattice fermions considered here, the critical polar- 
ization and the three-shell density profiles that we have 
found out are in good agreement with 3, ll| (and flilO . 
While the very large aspect ratio a ~ 45 used in [14| is 
not accessible because of computational limitations, we 
have not found any signatures of changing density profiles 
when increasing the trap aspect ratio. 

The exotic superfluid phases found in our study sug- 
gest that high polarizations P < P c can be an interesting 
area for future experiments with spin-imbalanced Fermi 
gases. In our calculations at T — 0, we have found that 
the FFLO-like phase emerges at high polarizations close 
to P c . The amplitude of the oscillating pair potential 
across the trap center is expected to be still significantly 
large, which implies that it may be experimentally ac- 
cessible at low but finite temperatures. In addition, the 



density and pair potential profiles that we have found in 
elongated traps reveal the 3D characteristics which are 
essentially different from the profiles predicted in strictly 
ID. This emphasizes the possibility for the experimental 
observation of the FFLO phase in 3D systems. 
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